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Abstract 

This paper addresses the problem of finding multiple near-optimal, spatially-dissimilar 
paths that can be considered as alternatives in the decision making process, for finding 
optimal corridors in which to construct a new road. We further consider combinations 
of techniques for reducing the costs associated with the computation and increasing 
the accuracy of the cost formulation. Numerical results for five algorithms to solve the 
dissimilar multipath problem show that a “bidirectional approach” yields the fastest 
running times and the most robust algorithm. Further modifications of the algorithms 
to reduce the running time were tested and it is shown that running time can be 
reduced by an average of 56% without compromising the quality of the results. 


1 Introduction 

In preliminary road design, selecting the best path for a new road is traditionally a long 
and political process. A variety of factors can require that road engineers design multiple 
alternative paths to be considered. Previous research has developed methods of modeling 
the road’s costs as well as computing the optimal alignment [l3l 02]. While this software 
provides useful insight for road engineers, it does not satisfy the need to hnd multiple road 
path options for review. There is a need for an algorithm that can efficiently compute 
multiple near-optimal, but distinctly different, paths that can serve as alternatives to the 
cheapest path. 
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This paper is concerned with the hrst step of the road design process: corridor selection. 
An initial path can then be rehned, optimizing the horizontal alignment within the corridor, 
using techniques from Mondal et ah [ST], Hare et ah or citations therein. 

To locate the initial corridor we model the terrain as a three dimensional grid of points to 
construct a spatial graph. Adapted forms of Dijkstra’s Algorithm are then used to select the 
set of spatially-dissimilar corridors, while minimizing an approximation of the earthwork and 
pavement costs. Additional techniques to reduce the running time have also been developed. 
These include an adapted form of the A* algorithm and two methods of imposing height 
restrictions, to avoid spending time searching for unrealistic roads that are far from the 
ground. 

1.1 Road Design 

Previous work on road alignment selection has developed a number of discrete models. The 
simplest and most common is a regular grid |13| , where the center of each grid cell is given a 
vertex and then edges are dehned by the vertices of the eight adjacent grid cells. This model 
restricts possible road trajectories to having eight directions. Another method, presented by 
Trietsch [121 E], used a honey-comb grid, which allows for angles in multiples of 30°. In 
many previous discrete models the space between each adjacent vertex was as small as 200 
m and up to 2 km [laimss], which does not allow for detailed and accurate assessment 
of costs in a given region of a cell. Many models also formulated the construction costs 
independently of the road’s direction [13] . 

On the horizontal alignment problem, Huber and Church ira took an in-depth look at 
minimizing the errors in the cost evaluation associated with path planning problems, such 
as road design, by increasing the number of possible directions of movement. Later, Lee et 
al. p3| applied a neighbourhood search technique to approximate horizontal alignments with 
a piecewise-linear curve. Once they selected a horizontal alignment they rehned it to be a 
smooth path meeting road standards such as curvature restraints. Easa and Mehmood m 
used collision frequency data for various road types to improve the safety of horizontal 
alignments. Kang et al. [2S] improved an existing Horizontal Alignment Optimization (HAO) 
model by restricting the search space with feasible gates. 

Kang izg used genetic algorithms to choose new highway segments that intersect with an 
existing highway network. A bi-level approach was developed to hrst select candidate paths 
with intermediate solutions to the genetic algorithm and then evaluate them for traffic how 
optimization [2^ \27\ . 

Jha [T9| approached the road design problem with genetic algorithms. Further work was 
done to integrate GIS to include the costs of land boundaries, environmental impact, to¬ 
pography, travel time, and noise and air pollution of highways near cities [21112Q|. Jha et 
al. [22] reviewed cost formulation as well as common road design optimization algorithms. 
Continuing this work Yang et al. im adapted the genetic algorithms to hnd alternate routes. 
Bosurgi et al. considered environmental constraints using particle swarm optimization [3]. 
They also added new types of curves and provided a genetic algorithm approach to optimize 
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the parameters [1]. Shafahi and Bagherian proposed a customized particle swarm opti¬ 
mization algorithm [30] while another 3D highway alignment model solved by evolutionary 
algorithms was considered by Jong and Schonfeld [2S]. It is also worth mentioning Jong’s 
Phd Thesis 

User interface designs have been developed and proposed for the road design problem. 
Church et ah [6] designed one using a multi-objective model to allow hnding both optimal 
and alternative paths. 

Many of the previous road design optimization methods are focused on hnding a single 
path or corridor. In practice, road design is a political process in which it is impossible 
to determine and evaluate all of the environmental and political cost factors ahead of time. 
Rather than using a multi-objective model with what will surely be an incomplete set of costs 
and constraints, this paper focuses on hnding k spatially dissimilar paths of nearly equal 
cost that can be reviewed by a board of engineers, politicians, and/or environmentalists for 
additional costs and impacts. 

1.2 /c-shortest-path Algorithms 

Extensive research has been done on the problem of hnding the fc-shortest-paths in a net¬ 
work [3 El HU Ea sg. In general there are two diherent approaches: deletion algo¬ 
rithms [U El EH] and deviation algorithms [ll |T6] . 

Deletion algorithms propose using a conventional path-hnding algorithm, such as Dijk- 
stra’s, to hnd the optimal path. The edges from the optimal path are successively deleted 
from the graph and the path-hnding algorithm is re-run to generate a set of secondary paths. 
The cheapest among them is selected and the process iterates. 

Deviation algorithms use the information generated by a shortest-path tree to the des¬ 
tination to exploit the frequent locality of the fc-shortest-paths m [16]. They begin with 
the cheapest path and then search for the deviation that ohers the smallest increase in cost. 
While this property of the /c-shortest-paths is what allows them to achieve the most compet¬ 
itive time complexity, it also offers insight as to why it, and indeed any of the conventional 
fc-shortest-path algorithms, are not well-suited for adaptation to the spatially-dissimilar mul¬ 
tipath problem. 

Both families of algorithms perform poorly when applied to the multipath road design 
problem, which corresponds to a fc-shortest-path problem applied to a dense spatial grapliQ 
with the additional complication of hnding spatially dissimilar paths. 

As an academic example consider the behaviour of a fc-shortest-path algorithm on a 
uniform cost grid (something that loosely approximates a very hat prairie). The globally 
optimal path would be the straight line pi seen in Figure [T} Paths p 2 , ps, and p 4 would all 
be equally priced and so any of them may be the second path found depending on how the 
algorithm breaks ties. Of the three options presented p 4 would rank as the most desirable 
since it is the most diherent from pi. However, p^ is still nearly identical to pi when using a 

^We call a dense spatial graph, a spatial graph with vertices having degree 26 except on the boundary. 


3 



refined grid and so none of these paths should be considered acceptable as candidate alternate 
paths to pi- 



Figure 1: A uniform cost grid and four possible paths illustrating the expected behaviour of 
a fc-shortest-path algorithm. The cheapest path is pi, and possible second, third, and fourth 
cheapest paths are p 2 , ps, and p^. 

The simplest way to adapt these algorithms to the multipath road-design problem is to 
iteratively generate paths until we have found k paths that satisfy our spatial dissimilarity 
criterion. Since early iterations of this approach would be spatially similar paths, this method 
would have to be iterated quite extensively before good alternative routes are found. As such, 
this approach performs poorly in practice. 

The dissimilar path problem is one often proposed in the context of transporting haz¬ 
ardous materials across highway networks PEsiian]. The goal of the dissimilar multipath 
problem is to hud a set of spatially dissimilar paths between a source and a destination. A 
variety of dissimilarity indices have been used to approach this problem. For two given paths 
many of them dehne a relationship between the shared length of the paths and the non-shared 
length of the paths PIMIEI]. This method is well-suited to the hazardous material (Haz- 
Mat) objective function which is based on the value of risk rather than the transportation 
costs. It yields local minima that are spatially diverse. Dissimilar path algorithms are better 
suited to hnding alternative paths for road design than standard fc-shortest-path algorithms. 
However, in the context of road design this method will hnd paths such as p^ seen in Figure 
These paths will share a minimal amount of length, but will be spatially very similar, due 
to the dense nature of the road design grid. Indeed, minimizing the shared length of paths is 
best suited to existing highway networks where alternate paths consist of different highway 
options. 

Dell’Olmo et ah [Sj, also in the context of HazMat transportation, introduces an inter¬ 
esting modihcation to previous dissimilarity indices by dehning a buffer of area around each 
of the paths and calculating an index determined by the area of overlap. This method is 
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particularly useful in HazMat transportation where the buffer area can be computed as the 
expected area that would be affected by an accident. 

Another metric used by Marti et ah |2Sl ESI ^ possible option that could have provided 
meaningful candidate alternative paths for road design. The distance metric used by Marti 
computes the average of the shortest distance from each vertex on path pi to p 2 and vice 
versa. The average distances of each path are normalized by the lengths of the respective 
paths and then the average of the two paths is taken as the dissimilarity index. 

While both of the previous two methods have the potential to provide meaningful candi¬ 
date alternative paths for road design, we opted to use another dissimilarity criterion, similar 
to one proposed by Lombard and Church [M]. Their metric looks at the area difference be¬ 
tween two paths using a method analogous to computing the absolute difference of the area 
under two curves. This method was chosen out of the three since it has equal potential to 
providing signihcantly different paths, while requiring the least computational expense. 

Most of the dissimilar path algorithms in the literature use a dissimilarity criterion related 
to the shared length of the paths and many only compare the alternate paths to the single 
cheapest path rather than to each other. This paper uses the area difference between two 
paths as a metric to ensure that the paths are spatially different. This type of metric 
is required for the dense graph used for the road design model, since it ensures that small 
deviations as seen in Figure are not accepted for alternate paths. This paper also measures 
the area between each path to ensure that no two alternate paths are spatially similar. 


1.3 Dissimilar Path Generation 


Several types of algorithms have been used to generate paths for the various dissimilar¬ 
ity criteria. A standard method among them is the Iterative Penalty Method [211231 EH]. 
Turner [13] applied a similar method to the multipath road design problem by increasing the 
weights of the edges of the original path. We adapt this idea for the road design problem 


with spatially dense graphs in Section [4^ 

Marti et ah [H51 ES] approach the problem using a Greedy Randomized Adaptive Search 
Procedure (GRASP) with Path Relinking. Jha [121 used Genetic Algorithms on the problem 


of road design and returned intermediate solutions as candidate alternate paths. Later, Kang 
et al. [271129] used Genetic Algorithms to generate path alternatives designing new highways 
that intersect with an existing highway network. Yang et ah m build on this work using a 
multi-objective model and modifying the Genetic Algorithm. Genetic Algorithms were also 
used by Zhang and Armstrong [32] for the multi-objective corridor problem. We did not study 
genetic algorithms, or GRASP, in this paper in favour of deterministic algorithms. Future 
research should explore comparisons between the approaches herein and nondeterministic 
methods, as well as explore the possibilities of hybrid approaches. 

A different method on multi-objective corridors was used by Dell’Olmo et al. [H] to hnd 
the pareto-optimal set of non-dominated paths. Despite the different objectives we suggest 


an algorithm of a similar nature in Section 4.3 


A common strategy is to generate a large set of candidate paths and then reduce this set 
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based on the chosen dissimilarity criteria. Kuby et ah [31] used a minimax method to select 
the set of paths which maximized the minimum dispersity between the paths, a problem 
which is known to be NP-hard ra- On HazMat transportation Kang et al. [201 generate 
their initial path set using regular fc-shortest-paths algorithm. This is appropriate for their 
problem since their objective function is based primarily on risk rather than travel costs. 
This objective function, unlike one based on construction costs, will often produce best 
paths that are spatially dissimilar. 

The Gateway Shortest-Paths method was introduced by Church and Lombard illSl] as 
a method of generating a set of spatially dissimilar candidate paths. We use a similar path 
generation method in Section 4.4 A very recent method was developed by Scaparra et al. [39] 
using a multi-gateway shortest-path method. This algorithm shows promise in generating 
a large number of candidate paths, but requires additional time, requiring a shortest-path 
algorithm to be run for every vertex in the graph. Given the density of vertices in the road 
design network this approach would take too long to run and would generate more paths 
than necessary to select reasonable alternatives. 

We introduce our road design model in Section and discuss two height restrictions to 
decrease running time in Sections 2A and 2.5[ Section [^contains the dissimilarity constraint 
required for our corridors. Five dissimilar multipath algorithms are presented and their 
theoretical performance is discussed in Section]^ Numerical results were collected and are 
summarized in Section Section contains some concluding remarks and ideas for further 
improvements in future work. 


2 The Model 

2.1 The 2T)-uj model 

In order to achieve accurate cost assessments and realistic paths we are using a very dense 
grid with approximately 10 meters between each vertex. However, with the increased detail 
of the model a new challenge arises: an alignment with a 90 ° turn (or more) between 
two vertices will have to take place in a very short space. This abrupt turn will generally 
violate road safety curvature constraints, making these alignments infeasible. We begin by 
introducing a constraint that the sharpest turn any road can take is a 45 ° angle. While this 
constraint may prove sufficient for some applications, e.g. logging roads, two consecutive 45 ° 
turns would still yield an unsafe highway alignment. However, this method is still sufficient 
for producing an initial corridor in which further curvature constraints can then be applied 
during the hne-tuning of the horizontal alignment optimization process. 

Simply storing the current direction and applying this restriction results in a constrained 
shortest-path problem. In order to use classic shortest path algorithms, such as Dijkstra’s 
Algorithm, we create a new 2D-a;’ model based on an augmented graph whose shortest-path 
solution corresponds to a solution to the constrained shortest-path problem. Trietsch [42] 
used this idea with hexagonal grids to apply curvature constraints to preliminary road de¬ 
signs. 
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Figure 2: The 8 possible ojh directions. 


The augmented graph has eight vertices for each {x,y) point, each one corresponding to 
an incoming edge orientation, see Figure A simple way of abstracting this information 
is to consider the direction of travel as another axis, let this be the cu/i-axis, for horizontal 
orientation. Each direction can then be assigned a coordinate, 0-7. 

Now that we have a new coordinate, we need to modify our original coordinate system. 
Let v{x,y,0Jh) be a vertex with Euclidean coordinates {x,y) and orientation Uh- A vertex 
with orientation Uh, will then have an outgoing edge to vertices with orientation o;^, where 

{ (^h - 1 , 

ujh, (mod 8), 

+ 1 , 

Theorem 2.1 Finding the shortest-path with the 2D-oj model can be solved in log-linear 
time. 

Proof: In the original model, every vertex has at most 8 edges, one to each adjacent vertex 
(with vertices on the boundary having slightly less edges). So if our map contains V vertices, 
we have at most 

E < 8E/2 = AV edges. 

In our new model, we have 8 different possible orientations per vertex. Let V be the new 
number of vertices, which is now 

V' = 8V. 

However, for this model, since we are restricting the possible directions of travel to 45° 
angles, we only have 3 possible outgoing and symmetrically 3 incoming edges at each of our 
augmented vertices, i.e., we have at most 

E' < 6V'/2 = 24V = 6E edges. 
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Thus, the number of vertices and edges is increased by a constant factor, which means that 
Dijkstra’s Algorithm will still run in log-linear time. ■ 

An added beneht of this model is that without any further modihcation, we can now 
easily improve the accuracy of our cost assessments. Since we are now storing the direction 
of travel, we can assign different costs to each edge, based not only on their location, but 
also their orientation. This increases the accuracy of the cost function, as road alignments 
traveling up a hill vary significantly in cost to those traveling along a hill. 

2.2 The 3D-a; model 

In order to extend our 2D-Ci; model to three dimensions, we add another axis, call it for 
vertical orientation. Unlike ooh, oOv only needs 3 orientations: -1, moving down; 0 maintaining 
elevation; and 1, moving up. Note v{x,y, z,Uh,uJv) the vertex with Euclidean coordinates 
{x, y, z) and orientation {uh, iXy). A vertex with orientation {uh, cZy) is then a vertex with an 
incoming edge with orientation (uh,cOy) and an outgoing edge to vertices with orientation 
where 


fcoh-1, 

u'he I Uh: (mod 8), 

[cu/i -|- 1, 

{ min(a;„ -h 1,1), 

UJv, 

max(a;^ — 1, —1), 

Theorem 2.2 Finding the shortest-path with the 3D-u model can he solved in log-linear 
time. 

Proof: Similarly to our 2D-a; model, the number of vertices and edges increases by a 

constant factor. The basic 3D model has at most 24 edges per vertex (with boundary 
vertices having slightly less). If our map contains V vertices, then we have at most 

E < 24U/2 = 12U edges. 

In the new model we have 24 different possible orientations per vertex, which increases 
our number of vertices to 

V' = 24V. 

Since we are restricting the possible directions of travel to 45° angles, we have up to 9 
outgoing and 9 incoming edges at each of our augmented vertices (vertices with = 0 have 
9 and vertices with = — 1 or 1 have 6). Which means we have at most 

E' < 18U72 = 9(24U) < 18E edges. 



The number of vertices and edges increases by a constant factor, which means that Dijkstra’s 
Algorithm will still run in log-linear time. ■ 


In practice Dijkstra’s Algorithm is able to hnd solutions, but requires signihcant com¬ 
putational time due to the size of the grid used. We discuss three methods of adapting 
Dijkstra’s Algorithm specifically for the road design problem in order to reduce the overall 
running time required. These modihcations can also be used in combination to speed up all 
of the dissimilar multipath algorithms in Section which use Dijskstra’s algorithm as their 
underlying shortest-path algorithm. 

The A* algorithm [7j is one that is commonly used as it maintains the same worst-case 
complexity as Dijkstra’s Algorithm, but in practice it has a faster expected running time. It 
uses a heuristic for a lower bound of the distance remaining to the destination and modihes 
the arc weights accordingly. For our application we are using the cost of paving a straight 
road to the destination. When used in the Bidirectional Selection Algorithm presented in 
Section 4.4 we modify the heuristic using the formula given by Ikeda et al. [IH], which 
provides the best known running time without losing global optimality. 


2.3 Cost Formulation 


In this model costs are assigned to each edge in the graph. Depending on the information 
available, a large variety of costs could be incorporated into the edge cost such as paving, 
earthwork, land acquisition, user costs, expected accident rates, etc. Some costs, such as 
expected accident rates or bridge construction, may be more complicated to quantify as the 
cost cannot be determined without considering the surrounding edges used in the road. 

For the numerical results in Section 5.2, we used only earthwork and paving costs. This 
approach was adopted to make edge costs easy to compute. The paving costs are simply 
proportional to the length of the edge, while the earthwork costs are based on the volume of 
earth excavation or embankment. Earth excavation and embankment was calculated using 
the height of the edge relative to the elevation profile of the ground. Namely, the ground is 
a piecewise linear function while the edge is linear so for each linear piece of the ground, we 
compute the area difference with the edge and multiply by the width of the road to obtain 
the volume. While changing cost structures may alter the hnal solutions determined, we 
believe that algorithmic performance would only be minimally effected. 


2.4 Simple Height Restriction (HR) 

We note that it is more expensive to build roads that are far above or below the ground than 
to build roads along the ground. This trend is not exploited by Dijkstra’s Algorithm which 
wastes time searching unrealistic paths. We introduce a simple novel constraint that restricts 
how far from the ground the algorithm can search. We note the following parameters 

- Hmax maximum displacement allowed from the ground, 

- R radius checked around a given point. 
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In the event the gronnd changes elevation abrnptly, we may need to build a slope to allow 
our road to change from one elevation to the other. As a result, for these parts of the map 
we will need to allow the algorithm to search farther away from the ground, which is why 
we include our radius R. Let Z{x,y) be the elevation of the ground at {x,y). Then for any 
position {x,y), a maximum distance above the ground, 6h+{x,y), is computed by 


where 


5h+{x,y) = max{max{Z{xi,yj)},Z{x,y) + H^nax) - Z{x,y), 

X — R<Xi<x + R 


Z{xi,yi) = <{xi,yj) : 


y-R<yj<y + R 

Symmetrically, we have the maximum distance below the ground 


5h_{x,y) = mm{m\n{Z{xi,yj)},Z{x,y) - Hmax) - Z{x,y). 


If R is set to 0, then only vertices within Hmax of the ground will be searched, which 
imposes the height restriction. The value of R increases the radius searched at each point 
when checking to see if the ground elevation abruptly changes. If the ground within the 
radius extends beyond Hmax, then the height restriction is relaxed to allow these points to 
be connected. 

As Hjnax and R decrease, significant reductions in computation time can be observed. 
However, one must be careful to not over-constrain the problem, which could remove desirable 
roads from being considered. 


2.5 Expanding Height Restriction (EH) 

The Expanding Height Restriction method was designed to initially only consider a small 
subset of vertices that are near the ground. A set of rules were devised to expand this subset 
when necessary to ensure that the problem was not over-constrained, removing promising 
candidate paths. 

In particular there are two cases where building a road that is not near the initial ground 
elevation may have been cheaper. The first is when a cliff steeper than the maximum 
permissible road grade needs to be climbed, which may require a “ramp”-like structure to 
build the road. The second corresponds to when it is cheaper to build a straight road through 
an obstacle rather than to incur the added pavement costs of building a longer road to avoid 
the obstacle. 

In numerical experiments (see Section]^, these rules had no effect on solution quality, 
but signihcantly slowed solution time. 


3 Dissimilarity Constraint 

We use an area constraint between each of the k alternate paths selected by the algorithms 
to ensure that the alternate paths are not small deviations of the globally optimal path. 
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Definition 3.1 We denote by SA{p,q) the percentage area difference of paths p and q. The 
percentage area difference between two paths is found by first projecting the paths onto an x-y 
plane, computing the area between the two paths, and then dividing by the area of a rectangle 
whose width is equal to the width of the map and whose length is equal to the straight-line 
distance between the endpoints of the paths. See Figure\^ 



Figure 3: The projection of the area between two paths onto an x-y plane, bounded by the 
endpoints of the two paths. The percentage of area is then calculated by dividing the area 
of the projection, by the area of the corridor between the two gray walls shown. 

We denote by SAmin the minimum percentage area difference required between each of 
the alternate paths selected. 

Lemma 3.2 The worst-case time complexity to compute the area between two paths is 0{L), 
where L is the number of vertices in the longest path found. Note that L n, where n is 
the number of vertices. 

Proof: The area computed between two paths is computed analogously to the integral 

of the absolute difference of the area under each path. This requires passing the length of 
the paths three times: once each to compute the area under each vertex and once more to 
compute the absolute difference. Since each pass requires 0{L) operations, the total time 
complexity is 0{L). ■ 

This dissimilarity constraint was selected both for its relatively low computation time and 
its simplicity in implementation. While there are many other possible choices for dissimilarity 
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metrics and this may not be the best one, this was believed to be an adequate metric to 
enforce spatial path dissimilarity. The constraint is related to the width of the map to help 
avoid over-constraining the problem. Should the dissimilarity constraint be set too high, 
most problems will become infeasible. By incorporating the width of the corridor into the 
constraint, future work can be done to determine a range of values which will typically yield 
feasible solutions with maximal dissimilarity. 


4 Dissimilar Multipath Algorithms 

This Section presents the dissimilarity criterion used by the hve dissimilar multipath algo¬ 
rithms presented in the following subsections. Pseudo code algorithms are given and their 
theoretical performance is discussed. 


4.1 Sensitive Elimination Method (SE) 

The Sensitive Elimination Method is an original algorithm similar in concept to deletion 
algorithms for Ending the fc-shortest paths [5l |32l US] . Instead of removing each of the edges 
in the optimal path, this method seeks to identify a promising edge to cut. The Sensitive 
Elimination Method begins by computing the optimal path using Dijkstra’s Algorithm. Next, 
some “sensitivity analysis” is performed for each edge in the optimal path, to determine 
which edge is the most sensitive. Then that part of the map is removed from the set of 
feasible edges and the shortest-path algorithm is again run on this reduced map. We use the 
following notation 

- w > 1 - sensitivity width, 

- k > 1 - target number of dissimilar paths, 

- 6Amin - minimum area difierence required from every other path as a percent, 

^ SCmax - maximum price difierence from the optimal path, as a percent. 

The sensitivity of a given incoming edge is computed by finding the difierence in price 
of building an edge with the same vertical position and orientation, {z,u)h,^v), but shifted 
to the left and right by d = 1, 2, ...w. This then produces an array of cost differences to the 
right and an array of cost differences to the left. Then, the sum of the absolute value of each 
array is multiplied to give a single numerical representation of the sensitivity of that edge. 

The cost model in this paper uses earthwork costs and paving costs. Since the segments 
to the left and right are the same length, the paving costs will be the same. The difierence 
in costs are based on the earthwork costs, which in turn is based on the displacement of the 
road segment from the ground. For simplicity, instead of computing the difierence in edge 
costs, we use a simple surrogate and compute the difierence in elevation at the vertex to 
which the edge points. 


Example 4.1 Consider computing the sensitivity with w = 2 of the edge going into vertex 
A. Let the black line in Figure fh represent the elevation profile to the left and right of 
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Figure 4: The black line is the elevation prohle of the ground to the left and right of vertex 
A. The solid gray line shows the height of the ground at A and the doted gray lines show 
the elevation differences to the left and right of A. 


vertex A. The elevation differences to the left and right are [—0.8,—0.5] and [0.1,—0.2] 
respectively. Taking the absolute value and summing the elements of the arrays gives 1.3 
and 0.3 respectively. Finally, we multiply these scores to get a raw index of 0.39, which 
can he compared with other sensitivity scores. We note that in this example the left side is 
considered more sensitive than the right, but the final sensitivity index is closer to that of the 
right. We chose to multiply the sensitivity values of the left and the right so that edges that 
are sensitive on both sides are selected first. 

For a given optimal path, let A be the vertex pointed to by the edge most sensitive to local 
perturbations. Instead of only removing A, we also eliminate all of the vertices with {x, y) 
coordinates equal to those of the edges to the left and right of A that were used to compute 
the sensitivity of A. By doing this, we have then created a wall in our map centered around 
the most sensitive vertices of our original path, through which any subsequent alternative 
paths cannot pass. See Figure]^ 

Algorithm: 

while we have not yet found k paths 
Compute the optimal path, 
if no path was found 

Replace the edges that were cut from the previous iteration. 

Select the next most sensitive, untried edge of the previous path found. 

Remove it and those to the left and right, 
elseif the new path is too expensive 

Stop, subsequent paths will also be too expensive, 
elseif the new path is too similar to one of the old paths 

Replace the edges that were cut from the previous iteration. 

Select the next most sensitive, untried edge of the previous path found. 
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Figure 5: Visualization of the Sensitive Elimination Algorithm. The optimal path’s most 
sensitive edge was chosen and the white wall shows the part of the map that was removed, 
making it impassable. The second path, on the left, is the resnlting alternate path fonnd. 


Remove it and those to the left and right. 

else 

Save the newly fonnd path. 

Compnte the sensitivity of each edge. 
Select the most sensitive edge. 

Remove it and those to the left and right 

end 

end 


One of the problems with this method, is that it may reqnire computing many paths that 
do not get used, a fact that is reflected in the worst-case time complexity. 

Proposition 4.2 The worst-case time complexity for the Sensitive Elimination Algorithm 
is 0{kLn\og{n) + where n is the number of vertices, k is the number of paths to be 

found and L is the number of vertices of the longest path considered, with L <^n. 

Proof: The algorithm Erst finds the optimal path containing np to L vertices. It is possible 
for all bnt one cnt to result in either an infeasible map, or an alternate path that is too similar 
to another path. This can then result in Dijkstra’s Algorithm being run np to L times to 
find a single feasible alternative path. Since we need to find k — 1 alternate paths, we may 
need to rnn Dijkstra’s Algorithm 0{L{k — 1)) times, assnming that only the iteration 
results in an acceptable path each time. 

The predicted best cnt is selected from each of the L choices, which reqnires L operations 
to select the maximnm valne. Dijkstra’s Algorithm needs to be rnn for each attempted 
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cut. Dijkstra’s Algorithm requires 0(nlog(n)), giving us 0{n\og{n) + L) to generate each 
possible alternate path. Each time a new path is found it needs to be compared with the 
0{k) other paths, requiring 0{kL) operations (see Lemma 3.2). Combining these gives us 
0{n\og{n) + {k + 1)L) operations each time a new path is considered. 

Since we may need to repeat Dijkstra’s Algorithm at most L(k — 1) times, this yields a 
total time complexity of 0{L{k — l)(?7,log(?7,) + [k + 1)L)), or more simply 0{kLn\og{n) + 

k^L^). m 


4.2 Iterative Penalty Adaptation (IPA) 

A common method for Ending alternative paths is to iteratively apply cumulative penalties 
along the shortest-paths [2l[23l[38]- First, a shortest-path algorithm is run to find the optimal 
solution. Weights along the shortest-path are then increased by a penalty and the shortest- 
path algorithm is run again to find a second alternative path. This process is repeated 
iteratively to find the desired number of paths. 

Due to the dense nature of our graphs if we were to apply this method directly we would 
find that the secondary path found would often simply be the original path shifted ten meters 
to the left or right to avoid the penalty, which does not provide us with sufficiently distinct 
alternate paths. We adapt existing methods by not only increasing the weight of the optimal 
path, but by also increasing the weight of a corridor-like buffer around the optimal path to 
discourage the alternative paths from being chosen too close to the original paths. This is 
achieved by adding a weight to the optimal path and decreasing that weight linearly on both 
sides of the path till it reaches zero. The effect of the buffer is to consider that any minor 
variation of the optimal path can be built by the engineer on the ground since our goal is to 
identify the corridor, not compute a precise path. Hence, the buffer forces any new solution 
to be clearly distinct. 

Iterative Penalty methods inherently have many variable options that will determine the 
performance of the algorithm. One could use either an additive or a multiplicative penalty. 
For simplicity we have opted to penalize all of the edges near the original paths with an 
additive penalty. We used a decaying penalty - that is, one with less weight farther away 
from the original paths. For simplicity we used a triangular decay shape with a fixed ratio 
of length to height and introduce the following notations 

~ Wp - the initial penalty width, 

~ Wmax - the algorithm terminates when Wp exceeds this value, 

- k > 1 - number of paths, 

~ SAjnin - minimum area difference required from every other path as a percent, 

^ SCmax - maximum price difference from the optimal path, as a percent. 

Each iteration of the algorithm involves running a shortest-path algorithm and applying 
the corridor weights. The new path found will not always satisfy the price and area con¬ 
straints. If this is the case we will either decrease or increase the penalty width, Wp. We 


15 



maintain a bracket beginning with a lower bound of zero, initially without an upper bound. 
When Wp is increased we will initially double its value until it needs to be decreased, provid¬ 
ing an upper bracket. Each time Wp is changed the distance between it and the appropriate 
bracket is halved and the new bracket values are updated accordingly. 

Algorithm: 

while we have not yet found k paths 
compute the optimal path, 
if the new path is too expensive 
try the middle of the Wp bracket 
if we have already tried the new Wp 
stop 

end 

elseif the new path is too similar 

if we have an upper bound for the Wp bracket 
try the middle of the Wp bracket 

else 

double the Wp 

end 

if we have already tried the new Wp, or Wp > Wmax 
stop 

end 

else 

add the new path to the set of k paths 

clear data about which values of Wp have been used 

end 

save that the current value of Wp has been used 

end 


Proposition 4.3 The Iterative Penalty Adaptation method has a worst-case running time 
of 0{k‘^n\og{n) log{Wmax) + k^Llog{Wmax))- 


Proof : Each shortest-path algorithm iterations takes 0{n\og{n)) operations. Each time 
a path is found it must be compared to the 0{k)) other paths, requiring 0{kL) operations 
(see Lemma 3.2). This gives a time complexity of 0{kn\og{n) -|- kfL) for each iteration. 

However, since each iteration may not result in a feasible path the algorithm may have 
to be run again. Initially the penalty width, Wp, does not have an upper bracket. As long as 
the path found is always too similar to the original paths it will continue to double Wp. Since 
Wp is doubled it will run at most 0(}og{wmax)) times before reaching Wmax and terminating, 
or decreasing and setting an upper bound. The upper bound will be bounded by Wmax and 
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since each time the size of the bracket is halved, the algorithm will hnd the next path, or 
run out of new Wp values, in at most 0(log(rymax)) iterations. 

Each time a new path is found the bracket values are cleared and the data about which 
Wp values have been tried before are reset. This means that the algorithm can spend up to 
0{log{wmax)) iterations again to hnd each of the k paths, requiring a total of 0(log{wmax)) 
iterations. 

Since each iteration requires 0{kn\og{n) + k‘^L), this yields a total time complexity of 
0(Pnlog(n) \og{wmax) + k^L\og{w„,ax))- ■ 

4.3 /c-shortest-paths Adaptation (KSPA) 

The /c-shortest-paths Adaptation |16] (KSPA) algorithm was designed to take advantage of 
work that has already been done by computing all k paths simultaneously, storing common 
information for each path only once. It is based on an adapted form of a fc-shortest-paths 
algorithm that can be derived from a generalization of Dijkstra’s Algorithm. 

When a suboptimal path is found to a vertex, instead of discarding it, we store both 
the information of the optimal path and the suboptimal path, until we have at most k 
paths to each vertex, one optimal, and k — 1 suboptimal paths. In order to maintain the 
desired properties of our alternate paths, we enforce a restriction that the area between each 
suboptimal path to the same vertex is at least SAmin- 

Once we have found our optimal path, we can then set the termination condition to be 
when we have found k = k paths to the destination, or when the cost of new paths is more 
than 6Cmax the price of the optimal path. We note the following parameters 

- K = k - number of paths found to each vertex, 

- SCmax - maximum price difference from the optimal path as a percent, 

- ^Amin - minimum area difference required between every path as a percent. 

Let 

- active heap be a heap containing vertices that are adjacent to the known shortest-path 
tree from the source, 

- verteXa be the source vertex, 

- vertexd be the destination vertex, 

- pathss-d be a collection of found paths from the source to the destination, 

- dirriy be the maximum y value (width of the corridor). 

Also let each vertex have k paths to it, each associated with its own cost. (Rather than 
each vertex storing only the single best path and cost.) 

Algorithm: 

push Vs onto the activcheap with a cost of 0 
while we have not yet reached the destination: 
current = cheapest vertex in activeheap 
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for each neighbour of current 

if the new path costs within bCmax of the cheapest path to neighbour 
if the new path is more than 5Amin from other paths to neighbour 
if we have not yet found k paths to neighbour 

push neighbour onto activcheap with its updated cost 
elseif this path is cheaper than the most expensive path to neighbour 
replace the most expensive path and push this one onto activCheap 

end 

elseif this path is similar to exactly one of the other paths 
if the new path is the cheaper of the two 

replace the old path and push the new one onto activcheap 
end 

end 

end 

end 

end 



Figure 6: A uniform cost grid illustrating the behaviour of a fc-shortest-path algorithm. 

One limitation of the /c-shortest-paths Adaptation algorithm is the behaviour that is best 
illustrated on near-uniform cost grids as may be found in some prairies. As an academic 
example, consider the hve paths pi, ... ps on a uniform cost grid, as seen in Figure]^ Let B 
be the parent of A on this path. The path is the cheapest path. Paths p 2 and p^ are too 
similar to ps and so will be rejected. Paths pi and ps are dissimilar from p^ and so would 
be acceptable alternate paths, however they are too similar to p 2 and p 4 respectively and so 
will also be rejected. When each path is very similar to the paths on either side of it, the 
algorithm is only ever to hnd one path to each vertex and hence will only hnd one possible 
road alignment. 
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Proposition 4.4 The k-Shortest-Paths Adaptation Algorithm achieves worst-case bounds 
of 0 {k‘^L n\og{n)), where n is the number of vertices, k is the number of paths found to each 
vertex, and L is the number of vertices of the longest path considered, with L <^n. 


Proof: Since we are now finding at most k paths to every vertex, instead of each vertex 
being used in at most 1 path, we can now have each vertex used in at most k paths. This is 
then similar to having nn vertices, which means that we will have 0{Kn\og{Kn)) iterations 
of the algorithm. 

In each iteration of the algorithm we need to compare the new path we have found with 
up to K other paths. From Lemma 3.2 we see that each iteration will then have an additional 
0{kL) operations. This gives us a total worst-case time complexity of log(«:n)), or 

more simply log(n)). ■ 


4.4 Bidirectional Selection Method (BDS) 

The Bidirectional Selection Method (BDS) is a simple extension of Bidirectional Dijkstra’s 
Algorithm, similar to the approach used by Lombard and Church |3l]. Instead of using the 
regular termination condition for the Bidirectional Dijkstra’s Algorithm, we simply continue 
to grow both ends until we have found k dissimilar paths. Alternatively, we terminate the 
algorithm if the cost of the new paths found exceeds the maximum cost restriction, since we 
will not hnd any cheaper paths beyond those. 

The current version of our algorithm uses a method of path selection that is relatively 
simple. Each time a new path is found it is compared to the set of accepted paths and 
is added to the set, replaces one of the paths in the set, or is rejected. A computational 
shortcut is used by not considering a path again after it has been rejected. While this 
shortcut may result in a desirable solution being missed, the additional computation time 
required to process all possible road alignment groupings would be infeasible. 

A similar method was proposed to hnd dissimilar paths in road networks mm- As they 
discussed, the s-d paths produced by this algorithm have the specihc property that given 
an arbitrary vertex A an s-d path is formed by concatenating the cheapest path from Vg to 
A and the cheapest path from A to Vd- The parameters are 

- k > 1 - number of paths, 

- dCmax - maximum price difference from the optimal path as a percent, 

- 6Amin - minimum area difference required from every other path as a percent. 

Let 

- Vg be the source vertex, 

- Vd be the destination vertex. 

Algorithm: 

push Vg and Vd onto their respective heaps with a cost of 0 
while there are still entries in the source or destination heaps 
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grow from whichever side has a smaller heap 
if the newly added vertex forms a new path 

if this new path is more than SAmin from the other paths selected 
if we have not yet selected k paths 

save the newly found path in our set of k paths 
elseif this path is cheaper than the most expensive path selected 
replace the most expensive path 

end 

elseif this path is similar to exactly one of the other paths 
if the new path is the cheapest of the two 

replace the old path with the newly found path, 
end 

end 

end 

if we have not yet found k paths and new paths are less than 5Cmax of the cheapest path 
push the neighbours of the newly added vertex onto the heap, as needed 

end 


Lemma 4.5 The Bidirectional Selection method computes at most n possible alternate paths 
with worst-case time complexity 0{n\og{n)). 

Proof: For a given vertex, we are computing the cheapest path to the source and the cheap¬ 
est path to the destination. That is, we compute the cheapest path from the source to the 
destination, that passes through each vertex. This is done by running Dijkstra’s Algorithm 
twice, once from the source and once from the destination, each taking 0{n\og{n)). By 
concatenating these results, we then generate one path for each vertex, taking an additional 
0{n). This gives a combined worst-case time complexity of 0(nlog(n) -|- n), or 0{n\og{n)). 


Proposition 4.6 The Bidirectional Selection Algorithm has worst-case bounds 0{n\og{n) + 
kLn), where n is the number of vertices, k is the number of paths, and L is the number of 
vertices of the longest path considered, with L <^n. 


Proof: From Lemma 4.5| our algorithm finds at most one different path per vertex, we 

are then selecting k paths from a set of at most n paths. Each new path is compared 
to our set of up to k accepted paths, which, from Lemma 3.2, requires 0{kL) operations 
for each new paths. With n new paths this gives us a total worst-case time complexity of 
0(nlog(n) -|- kLn). ■ 
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4.5 BDS-KSPA Hybrid Method 


This algorithm combines the two strategies by running the KSPA algorithm from both 
directions, as is done in the BDS method. We do make note that the number of paths 
found to a particular vertex from one direction, as used in the fc-shortest-paths Adaptation 
algorithm, need not be the same as the total number of paths found. The value of k can 
be thought of as the number of shortest-path trees made from both the source and the 
destination, while q controls the number of paths that are selected and returned from the 
paths generated by the shortest-path trees. A value of k = 1 reduces the algorithm to the 
Bidirectional Selection method, but a large value of k, will increase the running time of the 
algorithm. 

We omit the pseudo code algorithm for this section. The modifications provided in 


Section 4.3 do not affect the new termination and path selection process described in the 


pseudo code of 4.5 


Proposition 4.7 The BDS-KSPA hybrid method has 0{K^Lnlog{Kn) + ungL) as a worst 
case time complexity, where L is the number of vertices in the longest path considered, with 
L <t:n. 


Proof: Similarly to the KSPA algorithm each shortest-path tree generation requires 

0{K‘^Lnlog{n)) operations as seen in Proposition 4.4 


Since we are now hnding at most k, paths to each vertex we are generating a set of 0{Kn) 
paths. Since this algorithm uses the same path selection process as the BDS method each 


path will require 0{qL) operations (see Proposition 4.6). This gives a total of 0{nKqL) 
operations for the path selection process. 

Combining these results gives us a worst-case time complexity of 0{K‘^Ln\og{n) PnKqL) 
for the BDS-KSPA hybrid method. ■ 


5 Numerical Tests 

Numerical tests were run to compare algorithm quality in terms of time required and ability 
to hnd a valid solution. We look at which algorithm is most consistently able to hnd spatially 
dissimilar paths that are near the cost of the globally optimal path. The running time is also 
compared to select the fastest algorithm with the best results. The numerical tests are also 
used to measure the improvement in the running time gained with the height restrictions. 

A parameter required for each test is 6Amin, the difference threshold to state two given 
paths are ‘distinctly different’. This value should be selected based upon the number of paths 
desired to be found, k. We have chosen k = 3 and we are using 5Amin = 12%. Another 
parameter required for both is SCmax, which is the maximum cost difference allowed in 
secondary paths, as a percent of the optimal path. We have selected SCmax = 10%. 

In our numerical results an algorithm is considered to have found a solution if it meets 
three criteria. First, we require three paths with one of them being the optimal path; second. 
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the maximum path cost must be within the specihed range of the cheapest path; and third, 
the minimum percentage area requirement must be satished. 


5.1 Test Environment 


The algorithms perform differently based on the dimensions of the map (in terms of number 
of vertices). For both of our test sets we are using approximately 10 meters between each 
horizontally-adjacent vertex and 1 meter between each vertically-adjacent vertex. The hrst 
test set has three different map lengths of 40, 80, and 160 vertices, while the second test set 
has map lengths of 320 and 640 vertices. For each road length we have 7 maps with a length 
to width ratio of 2:1 and 3 maps with ratios 8:1, giving us a total of 30 maps in test set 1 and 
20 maps in test set 2. The initial terrain data was obtained from the USGS National Map 
Viewer [H]. The individual maps were found by sampling mountainous and prairie regions 
and selecting a variety of different vertical dimensions to ensure a good spread of map types. 
We provide a detailed breakdown of the terrain features present in each map in Appendix 
[A| The exact data used can be obtained by contacting the authors. 

All of the tests for the hrst test set were run on the Orcinus cluster, a 9600 core available 
through Westgrid |13]. Each test was run on a single core of an Intel Xeon X5650 six-core 
processor, running at 2.66GHz. The different map sizes were allotted different amounts of 
memory, up to the 24GB of RAM available per node. The second test set was run on Grex, 
a 3792 core cluster available through Westgrid. Again, each test was run on a single core of 
an Intel Xeon X5650 six-core processor, running at 2.66GHz. Each node had up to 96GB of 
RAM available. All of the code was written in MATLAB R2013b. Note that the algorithms 
are serial; the only reason we performed the tests on a cluster was to speedup running the 
algorithms on the full test set. 

We have selected hve algorithms for testing in this paper, and three modihcations that 
can be made to each to improve the performance. We include results for each of the basic 
algorithms without modihcation. We then selected the two best algorithms, and present the 
improvements gained by the three modihcations. First they are combined with the heuristic 
described in Section 2.2, adapted from the A* algorithm. Next they are combined with 


both the adapted form of the A* algorithm, and the Simple Height Restriction from Section 
2.4[ Further tests were done with both the heuristic from the A* algorithm, as well as the 


Expanding Height Restriction from Section 2.5 


When the algorithms are using the Simple Height Restriction we have set the following 
parameters: R = 3 and H^ax = 1- The Expanding Height algorithms will be using a value 
of Hi = 0.5. We used these parameter values as preliminary experimentation indicated they 
would provide good results. 

The Sensitive Elimination method’s performance is related to the choice of w, which 
determines how much of the map is removed at each iteration and 5Ajnin- If uj is too small, 
then it is unlikely that the alternative paths found will be at least 6Amin different. As a 
result, we have chosen to use a value that can be calculated hy w = round(5Amm • Wm — 0.5), 
where Wm is the width of the map. 
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Preliminary testing done for the Iterative Penalty Adaptation method achieved optimal 
performance with an initial penalty width value of Wp = 10%, which is the value we used for 
our tests. 

For the BDS-KSP hybrid method we used a value of k = 2, so that at most two paths 
would be found to each vertex and a value of g = 3 so that a total of 3 paths would be found. 

5.2 Results 

We summarize the results of our numerical tests with the two performance prohles in Fig¬ 
ures and A performance prohle [9] shows on the y-axis the portion of problems that an 
algorithm was able to solve, while the x-axis is the time it took to solve those problems as 
a ratio of the time it took the fastest algorithm to solve each problem. 


Comparison of Algorithms for Time 



Figure 7: The performance prohle of all the basic algorithms measured over time. 

We compare our algorithms using the total running time of each algorithm as the basic 
metric. We can see from our results in Figure that the Bidirectional Selection method was 
able to solve the most maps in the least time. The next best algorithm was the Iterative 
Penalty Adaptation algorithm, followed by the fc-shortest-paths Adaptation algorithm. The 
Sensitive Elimination method was often unable to hnd adequate solutions before reaching 
a predetermined time-out time. We discourage use of the latter three methods, as none of 
them have shown promising results within reasonable time frames. 
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In the first test set, the Bidirectional Selection and the hybrid methods were able to hnd 
solntions for 20 of the 30 problems. Of the remaining 10 problems, 9 were unsolved by any 
method, and 1 was solved by the Iterative Penalty Adaption method. Since the problems 
were designed using real-world terrain data and none of the algorithms found a solution, it 
is possible that the 9 unsolved problems do not have three feasible roads that satisfy the 
costing and percentage area difference constraints. As such, we have removed the 9 unsolved 
problems from our data analysis (including performance prohles in Figures]^ and [^. 


Height Restriction and A* Component for Time 



Figure 8: The performance prohle of the Bidirectional Selection algorithm and the Hybrid 
Method combined with the Simple Height Restriction and A* component, compared to just 
the A* component, and the base algorithms. 

While all of the algorithms were run again with the heuristic modihcations for efficiency, 
the hrst three continued to perform poorly so we do not include their results here. We found 
that in each case the height restrictions show improvements to the efficiency of the algorithms. 
For simplicity we omit the results of the Expanding Height Restriction in Figure]^ as they 
provided the same quality of solutions as the Simple Height Restriction, but the added 
overhead meant that the running time was not competitive. 
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Our fastest algorithm, the Bidirectional Selection method with the Simple Height Re¬ 
striction and the A* component, reported a solution to 20 of the 30 problems in the first test 
set. It was able to solve each of these 20 problems in just over five hours, with the smallest 
problem being solved in only twenty-three seconds. The same method required less than 4 
GB with an average of 20 minutes for the problems which had a map length of 40 vertices; 
those with map length of 80 vertices took 7 GB and an average of 1.5 hours; and the rest 
with under 19 GB and an average of 2.3 hours, on maps with a length of 160 vertices. 

The Bidirectional Selection method with the Simple Height Restriction and A* compo¬ 
nent was run again for the second test set. It was able to solve nine of the hrst ten problems 
with map length of 320 vertices, but only one of the last ten problems with 640 vertices. 
Maps with length of 320 vertices had an average running time of 15 hours. The largest 
problem with 320 by 160 by 195 vertices required the most memory at 82.2 GB. The longest 
running time of all of the problems was a little more than 28 hours. Of the problems with 
length 640 vertices, only the smallest one with 640 by 80 by 59 vertices was solved, using 
21.2 GB and just over 16 hours. The other nine problems with lengths of 640 vertices hit 
the memory limit of 92 GB. 

Note that these timings should not be viewed as absolute since coding the algorithms in 
G instead of MATLAB will greatly reduce the running times. These results show the relative 
performance of each algorithm, in order to select the best one for further rehnement. 

5.3 Example Solution 

In this section we present a sample solution found by the Bidirectional Selection Method 
on a map with 320 by 160 vertices (3.2 km by 1.6 km), see Figure]^ The globally optimal 
path, pi, was determined to cost approximately 180,000 monetary units. The two alternative 
paths, p 2 and pa cost 190,000 and 193,000 monetary units, respectively. These costs satisfy 
the maximum cost constraint as they are both within 110% of the cost of the globally optimal 
path, at 105% and 107% respectively. The projected area difference between pi and p 2 is 
15%, the area between p 2 and ps is 13%, and the area between pi and ps is 28%. This meets 
the minimum area difference constraint which requires each path to be separated by at least 
12% of the area of the map. 
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Figure 9: The three paths returned by the Bidirectional Selection Method on a map with 
320 by 160 vertices (3.2 km by 1.6 km). This terrain corresponds to Test 9 in the second 
test set. 


6 Conclusion 


The primary contributions of this paper are the adaptations of several path-hnding algo¬ 
rithms to develop dissimilar multipath algorithms. The Sensitive Elimination algorithm in 


Section 4.1 is an original method belonging to the deletion family of /c-shortest path al¬ 
gorithms. The Iterative Penalty Adaptation algorithm in Section 4.2| modihes a common 
technique by adding a corridor to which a penalty is applied. All of the algorithms were 
modihed to include an acceptance test based on the dissimilarity criterion in Section We 
introduce a possible hybridization of the two methods BDS and KSPA in Section |4.5[ 

Our numerical results indicate that, of the algorithms tested, the only reasonable algo¬ 
rithm for practical use is the Bidirectional Selection method using A* and the Simple Height 
Restriction. The Bidirectional Selection method and its hybrid with the fc-shortest-paths 
Adaptation algorithm, were observed to outperform the other methods in both time and 
quality of solutions. The main restriction of the method is that it only considers paths that 
are a concatenation of the cheapest path from the source to a node and then the cheapest 
path from that node, to the destination. A result of this, is that its paths cannot begin by 
overlapping, then split, merge and split again. This property does not, however, stop the 
algorithm from producing two candidate paths which cross at an angle. We did not End this 
property to be a signihcant problem for our application, since it was still able to generate a 
large number of candidate paths. In practice, if an algorithm with this path property was 
unsatisfactory, then a hybrid of the Bidirectional Selection Method and the fc-shortest-paths 
Adaption algorithm could be used as described in Section |4.5[ While this does theoretically 


increase the range of paths that can be found our numerical results showed little difference 
in the quality of the solutions between the two algorithms. 
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Current implementations of the Bidirectional Selection method do not consider a path 
once it has been rejected. Future work may be able to improve the quality of the solutions by 
relaxing this condition. Such modihcations may also observe performance cuts and memory 
challenges that arise from storing and processing the additional paths’ information. Other 
possible improvements to the quality of solutions could be made by examining the types of 
solutions returned by the algorithms when using different types of dissimilarity metrics and 
constraints. With several other types of dissimilarity constraints available, it is possible that 
these metrics may yield substantially different set of solutions. 

The A* component of the algorithms combined with both the Simple Height Restriction 
and the Expanding Height Restriction reduced the computational time of the algorithms. 
The results from the Expanding Height Restriction are not included in this paper as they 
were not competitive with the Simple Height Restriction. Both height restrictions were able 
to solve the same problems, but the additional overhead required for the Expanding Height 
Restriction meant that it took longer. In some cases the height restriction algorithms were 
able to hnd more solutions within the preset time-out time than their regular counterparts. 
Our aggressive choice of height restriction parameters for the Simple Height Restriction pro¬ 
vided a signihcant improvement in time without decreasing the quality of solutions. Another 
appeal of the Simple Height Restriction is its simplicity in implementation. 

Before the Bidirectional Selection Method could be efficiently used in practice further 
rehnements would be necessary. It is worth noting that the largest test problem we considered 
was only 3200 meters. The time to solve such problems was approximately several hours. 
One possible direction for improvement would be to massively parallelize the algorithm to 
reduce the running time. 

Improving the accuracy of the heuristic used for the A* component of the algorithm 
would decrease the overall running time. The heuristic could be improved by the addition 
of unavoidable cost factors such as rivers, or mountain chains separating the source and the 
destination. 

Future versions of the dissimilar multipath algorithms could be modihed to connect two 
sections of an existing road network. This process would be straight-forward assuming the 
cost of connecting the new road to a specihc point of an existing network is known. Instead 
of simply adding one starting vertex, all candidate vertices along the existing road could 
be added with their respective costs. The destination nodes code similarly be modihed to 
trigger the stopping condition. Existing road networks bisecting the search space could also 
be incorporating by modifying the terrain graph to include low-cost or no-cost edges, as 
appropriate. 
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Appendices 

A Test Sets 

We have defined three terrain types to categorize the test maps to ensure that we have a 
diverse set of tests. With the interests of road construction in mind we have based this upon 
our maximum road grade of 10%. Let rrix^y be the maximum grade of the terrain in each of 
the eight Uh directions allowed for edges at the point {x,y). Then for each position {xi,yi) 
corresponding to the location of a vertex in our map, we calculate 

Definition A.l Let T{xi,yi) be the type of terrain at position {xi,yi). 

(a, if mxi,yi < 10%, 

T{xi, yi) = I B, if 10% < mx^,y^ < 20%, (1) 

[c*, if 20% <mx„yi. 

Using these definitions we then calculated the percentage of each category present in each 
test map. 
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A.l 


Test Set One 



Dim X 

Dim y 

Dim z 

A 

B 

C 

Test 1 

40 

20 

5 

100% 

0% 

0% 

Test 2 

40 

5 

9 

93% 

7% 

0% 

Test 3 

40 

20 

23 

54% 

27% 

19% 

Test 4 

40 

20 

37 

12% 

56% 

32% 

Test 5 

40 

5 

38 

0% 

0% 

100% 

Test 6 

40 

20 

53 

6% 

43% 

52% 

Test 7 

40 

20 

60 

0% 

8% 

92% 

Test 8 

40 

5 

60 

15% 

15% 

70% 

Test 9 

40 

20 

88 

0% 

20% 

80% 

Test 10 

40 

20 

90 

0% 

1% 

99% 

Test 11 

80 

40 

11 

91% 

9% 

0% 

Test 12 

80 

10 

10 

92% 

6% 

2% 

Test 13 

80 

40 

28 

37% 

31% 

32% 

Test 14 

80 

40 

39 

46% 

43% 

11% 

Test 15 

80 

10 

41 

29% 

56% 

16% 

Test 16 

80 

40 

67 

5% 

44% 

52% 

Test 17 

80 

40 

84 

38% 

16% 

47% 

Test 18 

80 

10 

61 

18% 

22% 

60% 

Test 19 

80 

40 

101 

17% 

21% 

61% 

Test 20 

80 

40 

123 

54% 

4% 

43% 

Test 21 

160 

80 

5 

100% 

0% 

0% 

Test 22 

160 

20 

14 

100% 

0% 

0% 

Test 23 

160 

80 

36 

71% 

19% 

10% 

Test 24 

160 

80 

50 

69% 

22% 

9% 

Test 25 

160 

20 

31 

69% 

28% 

3% 

Test 26 

160 

80 

85 

49% 

30% 

21% 

Test 27 

160 

80 

97 

41% 

30% 

30% 

Test 28 

160 

20 

27 

87% 

13% 

0% 

Test 29 

160 

80 

127 

46% 

19% 

35% 

Test 30 

160 

80 

146 

16% 

37% 

46% 


Table 1; The dimensions and type of terrain for the thirty tests used in the hrst test set. 
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A.2 Test Set Two 



Dim X 

Dim y 

Dim z 

A 

B 

C 

Test 1 

320 

160 

9 

99% 

1% 

0% 

Test 2 

320 

40 

20 

100% 

0% 

0% 

Test 3 

320 

160 

41 

75% 

21% 

4% 

Test 4 

320 

160 

67 

67% 

20% 

12% 

Test 5 

320 

40 

47 

78% 

10% 

12% 

Test 6 

320 

160 

111 

49% 

24% 

27% 

Test 7 

320 

160 

122 

39% 

34% 

28% 

Test 8 

320 

40 

132 

55% 

17% 

28% 

Test 9 

320 

160 

165 

39% 

24% 

37% 

Test 10 

320 

160 

195 

41% 

23% 

36% 

Test 11 

640 

320 

49 

71% 

22% 

6% 

Test 12 

640 

80 

44 

94% 

5% 

1% 

Test 13 

640 

320 

107 

68% 

18% 

14% 

Test 14 

640 

320 

192 

63% 

19% 

18% 

Test 15 

640 

80 

188 

64% 

24% 

12% 

Test 16 

640 

320 

274 

32% 

21% 

47% 

Test 17 

640 

320 

319 

20% 

27% 

53% 

Test 18 

640 

80 

332 

6% 

27% 

66% 

Test 19 

640 

320 

416 

7% 

17% 

75% 

Test 20 

640 

320 

479 

18% 

29% 

53% 


Table 2; The dimensions and type of terrain for the twenty tests used in the second test set. 
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